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1.0  SUMMARY 


This  report  describes  the  application  of  a  hybrid  large-eddy  simulation  /  Reynolds- 
averaged  Navier-Stokes  (LES/RANS)  turbulence  modeling  strategy  to  predict  reactant  mixing 
and  shock-train  propagation  within  a  model  Chemical  Oxygen  Iodine  Laser  (COIL)  unit.  The 
configuration  consists  of  a  converging-diverging  nozzle,  a  lasing  cavity,  and  a  diffuser.  The 
intent  of  this  investigation  is  to  assess  the  potential  of  LES/RANS  technology  in  predicting  the 
possible  effects  of  shock-train  motion  on  the  flow  within  the  lasing  cavity.  Results  have  been 
obtained  for  several  grids  and  for  several  back  pressures,  the  variation  of  which  fixes  the  average 
shock-train  position  within  the  nozzle.  Predictions  of  wall  pressure  are  in  reasonable  accord 
with  experimental  observations  for  both  LES/RANS  and  RANS  model  but  tend  to  under-predict 
the  initial  rate  of  pressure  increase  as  the  flow  encounters  the  leading  edge  of  the  shock  train. 
Though  the  shock  train  induces  a  transition  to  turbulence  as  well  as  local  flow  separation,  the 
effects  of  resolved  turbulence  and  unsteady  separation  do  not  appear  to  influence  the  mixing 
process  in  the  lasing  cavity  significantly.  Moderately  unsteady  laminar-flow  mixing  processes 
dominate  in  the  lasing  cavity,  providing  efficient  mixing  of  reactants  without  turbulence.  This 
indicates  that  current  diffuser  designs  may  provide  sufficient  isolation  of  the  lasing  process  from 
the  pressure  recovery  process,  at  least  for  moderate  back  pressures. 

2.0  INTRODUCTION 

Chemical  Oxygen  Iodine  Laser  (COIL)  systems  require  low  pressures  for  effective 
operation,  and  to  recover  ambient  pressure  levels,  diffuser  and  ejector  technologies  are 
necessary.  The  design  of  an  effective  diffuser  system  capable  of  containing  the  ‘shock  train’ 
that  provides  initial  pressure  recovery  involves  several  considerations.  Among  these  are  the  need 
to  operate  effectively  over  a  range  of  back  pressures  and  the  requirement  that  the  lasing  cavity  be 
isolated  from  the  unsteady  dynamics  of  the  shock  train  to  improve  mixture  homogeneity  and 
thereby  enhance  laser  gain.  Though  nearly  laminar-flow  conditions  may  be  maintained  in  the 
lasing  cavity,  a  transition  to  turbulence,  enhanced  by  shock-induced  amplification  of  boundary- 
layer  disturbances,  will  take  place  in  the  pressure-recovery  section.  Local  fluctuation  intensities, 
scaling  with  the  shear  layer  thickness,  emerge  as  a  result  of  turbulent  flow,  but  lower-frequency 
oscillations  also  can  emerge  and  can  lead  to  large-scale  oscillations  of  the  entire  shock-train 
system.  The  cause  of  such  oscillations  is  a  matter  of  debate,  even  for  single  shock  /  boundary 
layer  interactions  [1],  but  the  prevailing  opinion  is  that  they  result  from  net  entrainment  of  fluid 
into  separation  regions  at  a  frequency  generally  proportional  to  the  mean  size  of  the  separated- 
flow  region.  Gieseking  and  Edwards  [2]  show  that  the  dominant  low  frequency  component  of 
separation-shock  unsteadiness  is  directly  correlated  with  the  most  probable  residence  time  of 
fluid  within  the  separated-flow  region  -  an  observation  that  may  make  it  possible  to  predict 
unsteady  characteristics  by  extracting  residence  time  distributions  from  steady  calculations.  In 
COIL  laser  systems,  low-momentum  regions  are  found  near  comers  as  well  as  in  regions  directly 
affected  by  shock  impingement,  and  large-scale  oscillations  of  the  diffuser  shock-train  system 
are  indeed  observed.  It  is  of  interest  to  understand  the  causes  behind  shock-train  unsteadiness  in 
these  systems  and  their  concomitant  influence  on  the  lasing  process. 

The  objective  of  this  research  effort  is  to  apply  a  recently-developed  hybrid  large-eddy 
simulation  /  Reynolds-averaged  Navier-Stokes  (LES/RANS)  turbulence  model  to  simulate 
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shock- train  formation  in  a  model  COIL  nozzle  /  diffuser.  The  LES/RANS  method,  enhanced  by 
low-dissipation  numerics,  should  be  eapable  of  eapturing  the  larger  seales  of  turbulent 
unsteadiness  resulting  from  shock  /  viscous  layer  interactions  as  well  as  oseillations  of  the  entire 
shoek-train  system.  Predietions  are  compared  with  available  experimental  data  [3],  which 
consists  of  wall  pressure  distributions  obtained  in  a  representative  COIL  unit  eonsisting  of  a 
converging  /  diverging  nozzle,  several  injeetors,  a  lasing  cavity,  and  a  diffuser.  Previous 
predietions  using  GASP  on  a  33  million  cell  mesh  that  resolves  several  injection  units  have 
shown  a  signifieant  sensitivity  of  the  diffuser  response  to  variations  in  the  chosen  turbulenee 
model  for  a  given  baek  pressure.  [3,4]  The  remainder  of  this  report  describes  some  details  of  the 
numerieal  methods  used  to  sustain  turbulence,  the  LES/RANS  closure  model  employed  in  this 
work,  and  results  from  several  studies  that  foeus  on  idealized  and  realistic  COIL  nozzle  /  diffuser 
systems. 

3.0  METHODS,  ASSUMPTIONS,  AND  PROCEDURES 

3.1  Time  Advancement 

North  Carolina  State  University’s  (NCSU)  REACTMB  flow  solver  is  used  in  the  present  effort. 
REACTMB  solves  the  Navier-Stokes  equations  governing  a  multi-component  mixture  of  gases 
on  simply-eonnected,  multi-bloek  struetured  meshes  using  finite-volume  methods.  REACTMB 
is  designed  for  use  on  massively-parallel  computers  and  uses  MPI  for  message-passing.  Large 
struetured  meshes  are  deeomposed  into  a  number  of  smaller  blocks,  which  are  then  partitioned 
over  the  number  of  requested  cores  using  a  simple  divide-and-conquer  strategy.  REACTMB 
discretizes  the  Navier-Stokes  system  in  time  using  a  Crank-Nicholson  approaeh: 

Tjn+\,k  _Tjn  -I  1 

+-(1  +  +  -  (1  -  )  =  0  (1) 

where  Q  is  the  eell  volume.  At  is  the  time  step.  U  is  the  vector  of  conserved  variables,  and  R  is 
the  residual  vector.  The  funetion  6  is  defined  as 

e  =  -{\-  tanh(^^^)),  d,=\x  lO'V,  Ad  =  0.2d,  (2) 

2  Ad 

Here,  d  is  the  distanee  to  the  nearest  solid  surfaee.  The  function  6  switches  the  time 
diseretization  from  Crank-Nieholson  to  Euler  implicit  for  mesh  cells  essentially  within  the 
laminar  sub-layer.  Some  loss  of  temporal  aceuraey  results,  but  this  approach  is  necessary  to 
suppress  oscillations  in  the  pressure  and  transverse-velocity  fields  for  mesh  eells  with  a  very  high 
aspeet  ratio.  The  matrix  system  that  results  from  linearizing  (1)  is  solved  approximately  using 
a  bloek  incomplete  lower-upper  (ILU)  factorization  method,  and  the  system  is  converged  to  a 
prescribed  tolerance  over  a  sequence  of  sub-iterations.  Jaeobian  matrix  elements  are  stored  over 
the  number  of  bloeks  mapped  to  a  partieular  processor,  allowing  the  “freezing”  of  the  matrix 
elements  and  their  faetorization  over  the  duration  of  the  sub-iterations.  This  reduees  the 
computational  workload  significantly. 

3.2  Flux  Formulation 

To  sustain  turbulence,  it  is  necessary  to  reduce  numerieal  dissipation  signifieantly.  The  strategy 
employed  in  REACTMB  eombines  a  variant  of  the  Pieeewise  Parabolie  Method  [5]  with  a 
fourth-order  eentral  difference  seheme.  Edwards’s  low  diffusion  flux-splitting  scheme  (LDFSS) 
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[6]  is  used  as  the  Riemann  solver.  The  primitive-variable  vector  W  =  [p^,u,v,w,T,k,(X>Y  is  used 
in  the  reconstruction.  The  initial  step  in  the  PPM  reconstruction  sets  left-and  right  states  to 

7  .  1 


^LJ+l/2  ^R,i+l/2 


(3) 


which  yields  a  fourth-order  central  difference  approximation  on  uniform  meshes.  Monotonicity 
is  enforced  by  a  cell-by-cell  resetting  of  left  and  right  states  [5],  leading  to  final  left-  and  right- 
state  values  (superscript  ‘M’)  that  are  different  from  the  averaged  ones  (superscript  ‘A’).  The 
amount  of  numerical  dissipation  added  at  a  cell  interface  is  proportional  to  the  difference  in  left- 
and  right-state  values.  While  enforcing  monotonicity  helps  in  shock-capturing,  it  compromises 
the  ability  of  the  scheme  to  resolve  small-scale  turbulent  structures.  One  means  of  alleviating 
this  problem  is  to  blend  the  averaged  values  (3)  with  the  monotonicity-preserving  values,  so  that 
the  former  is  used  in  regions  of  high  vorticity  (boundary  layers,  shear  layers)  and  the  latter  is 
used  in  more  ‘inviscid’  regions,  where  strong  shocks  might  be  present.  A  function  due  to 
Ducros,  et  al.  [7],  defined  at  a  mesh  cell  as 

(v-vy 


f  = 


(V-F)"+| 


(0 


£•  =  1x10  /  max(Ax,  Ay,  Az) 


(4) 


is  used  for  this  purpose.  Here,  the  divergence  of  velocity  is  compared  with  the  vorticity  value.  If 
the  latter  is  much  larger  (in  shear  and  boundary  layers,  for  example),  the  function  moves  toward 
zero,  and  in  free-stream  regions  near  shocks,  the  function  approaches  one.  At  a  particular  cell 
interface  i+1/2,  we  use  the  function  as  follows: 


^L,i+l/2  ~  ^L,i+l/2  +  max(/;.,/;.^i 

^RJ+l/2  ~  ^R,i+l/2  +  max(y;.,/;.,„0.1)(IF, 


,0.mL%n  ^L,M/2^ 


M 

R,MI2 


■w, 


R,i+l/2 


) 


(5) 


This  scheme,  denoted  as  LD-PPM  for  low-dissipation  PPM,  is  used  for  all  LES/RANS 
calculations  presented  in  this  paper.  Viscous  and  diffusive  terms  appearing  in  the  equation 
system  are  discretized  using  second-order  central  differences. 


3.3  LES/RANS  Model 


NCSU’s  hybrid  LES/RANS  methodology  [2,8]  is  used  as  the  baseline  for  the  current  work.  In 
the  LES/RANS  model,  the  effects  of  anisotropic  near-wall  eddies  are  modeled  using  RANS 
concepts  (Menter’s A:-® baseline  model),  whereas  the  larger  turbulent  eddies  away  from  solid 
surfaces  are  captured  using  a  large-eddy  simulation  method.  The  shift  between  the  closure 
models  is  facilitated  by  modifying  the  eddy  viscosity  field  according  to 


k 

1 

f  .  1  1 

II 

(i-rK.,  +r- 

CO 

l-tanh[15(— -1)] 

1  '^N  y 

Here,  T  is  a  time-dependent  blending  function  that  connects  the  RANS  and  LES  branches  and  the 
quantity  2^  is  a  ratio  of  outer-  to  inner-layer  turbulence  length  scales,  with  the  former  calculated 

using  both  ensemble-averaged  and  instantaneous  turbulence  data.  Specific  definitions  are  as 
follows: 

=  giLter  Pouter  / hnner  ’  giLter  )  =  min[E)i ,  max(l,  ^  J^)]  (7) 

Lter  V 

where 
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/  =  c 

^  outer  ^  N  - 


1  Ovco  k  kn 


C^'^coco 


/.  =  Kid 

inner 


The  model  eonstants  Z)j  and  are  assigned  values  of  10  and  0.5,  based  on  calibrations  for  flat- 
plate  boundary  layers.  The  mesh  scale  is  taken  to  be  the  maximum  spacing  over  all  three 

coordinate  directions.  This  form  serves  to  shift  the  closure  to  unsteady  RANS  when  there  is  no 
possibility  of  resolving  the  largest  turbulence  length  scales.  If  the  maximum  mesh  scale  is 
selected  as  the  outer-layer  scale  consistently,  then  the  LES/RANS  model  behaves  similarly  to 
detached-eddy  simulation,  serving  primarily  to  isolate  turbulent  boundary  layers  from  massively- 
separated  regions.  The  quantity  is  the  resolved  turbulence  kinetic  energy,  calculated  via 
ensemble-averaging  as 

P^R  =  \ (8) 

Ip 

The  quantities  k  and  W  are  ensemble-averages  of  the  modeled  turbulence  kinetic  energy  and 
specific  dissipation  rate  variables,  which  are  obtained  from  Menter’s  model.  The  use  of 
ensemble-averaged  as  well  as  instantaneous  data  allows  the  RANS-to-LES  transition  location  to 
‘float’  about  a  time-mean  location  that  is  a  function  of  the  ensemble-averaged  state  of  the 
boundary  layer.  This  enables  the  blending  function  to  respond  more  directly  to  large  changes  in 
the  turbulence  length  scales  that  can  result  from  shock  interactions.  Ensemble-averages  are 
currently  calculated  using  an  exponentially-weighted  moving  average:  Q"  =  (1  - A)Q +  AQ" 
with  A  =  lAtl  T .  The  time  scale  r  is  defined  as  follows  for  the  cases  considered  herein: 


where  =  Llu^  is  defined  in  terms  of  the  length  of  the  domain  L  and  the  free-stream  velocity 
This  form  assumes  that  a  statistically-stationary  state  will  emerge  after  about  four  residence 
times 

4.0  RESULTS  AND  DISCUSSION 

This  model  has  been  applied  to  three  nozzle  /  diffuser  geometries  representative  of  COIL  lasing 
units.  The  first  contains  ~75  million  cells  but  neglects  injector  units  and  contains  some 
geometric  simplifications.  The  second  contains  ~70  million  cells,  incorporates  several  iodine 
injector  units,  and  corresponds  geometrically  to  an  experiment  described  in  [3].  The  third  is 
equivalent  to  the  second  but  with  many  more  cells  (228  million  cells  in  total)  added  to  resolve 
eddy  dynamics  better. 

4.1  Results  and  Discussion:  1®*  geometry 

Wall  pressure  distributions  for  two  back  pressures  (14.7  Torr  and  26.4  Torr)  for  the  first 
geometry  are  shown  in  Figure  1,  as  are  wall  pressure  standard  deviation  distributions,  a  measure 
of  shock-train  unsteadiness. 
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Figure  1:  Mean  and  rms  wall  pressure:  Geometry  #1 


Agreement  with  experiment  for  mean  pressures  is  fair  due  to  the  geometrie  simplifieations 
employed  in  the  first  geometry,  while  the  wall-pressure  standard-deviation  distributions  show  the 
impaet  of  shoek-train  unsteadiness  on  the  flow  in  the  lasing  cavity  as  well  as  in  the  diffuser 
section. 
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Figure  2:  Centerplane  Mach  number  distributions  at  different  time 
instances:  14.7  Torr  back  pressure 
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Figure  2  shows  centerplane  snapshots  of  the  shock  train  and  its  motion  for  the  14.7  Torr  back 
pressure.  Animations  of  the  pressure  signal  and  the  ‘footprints’  of  wall  separation  regions  reveal 
that  regions  of  reversed  flow,  located  primarily  in  the  comer  and  sidewall  regions,  enlarge  and 
shrink  in  response  to  shock  train  motion,  providing  a  path  for  disturbances  to  propagate  upstream 
to  the  lasing  cavity. 


Figure  3:  Centerplane  Mach  number  contours:  Geometry  #2 


4.2  Results  and  Discussion:  2"**  geometry 

Results  from  simulations  of  shock-train  development  within  the  second  geometry  at  a  back 
pressure  of  22.8  Torr  are  shown  in  Figures  3-7.  Figure  3  shows  centerplane  Mach  number 
contours  for  two  cases:  a  steady  RANS  simulation  and  a  time-dependent  LES/RANS  simulation. 
Some  eddy  content  in  the  diffuser  is  noted,  but  the  stractures  appear  to  be  too  large,  indicating 
that  the  mesh  is  too  coarse  in  this  region. 


{Turt.  Length Vc()i:  0  0.2  0,4  0.6  O.a  1  1,21.41.01,0  2 


Figure  4:  X-Y  centerplane  length-scale  ratio  from  RANS  solution:  Geometry  #2 


Figure  4  plots  the  ratio  of  the  modeled  turbulence  length  scale  from  the  RANS  solution  to  the 


maximum  of  the  grid  dimension:  ratio  = 


Adequate  resolution  would  imply 

max(Ax,  Ay,  Az) 

that  this  quantity  is  greater  than  one  in  regions  of  significant  turbulence  activity.  Figure  5 
compares  this  ratio  along  the  centerline  of  the  X-Y  plane  with  that  extracted  from  the  Geometry 
#1  RANS  solution.  In  the  diffuser  section,  the  ratio  is  always  greater  than  one  for  Geometry  #1 
while  generally  less  than  one  for  Geometry  #2.  This  indicates  that  the  LES/RANS  calculation 
for  Geometry  #2  is  under-resolved  in  this  region  and  in  other  locations.  As  a  possible 
consequence  and  as  shown  in  Figure  6,  the  wall  pressure  prediction  is  very  poor  for  this  case. 
Figure  7  shows  some  evidence  of  unsteadiness  of  the  mixing  plume  in  the  lasing  cavity.  This 
may  arise  from  natural  instabilities  that  possibly  are  amplified  by  external  forcing  generated 
from  shock-train  motion. 


Approved  for  public  release;  distribution  unlimited. 

6 


Figure  5:  X-Y  centerline  length-scale  ratios: 
Geometry  #1  vs  Geometry  #2 


Figure  6:  Wall  pressure  distribution  (22.8  Torr  back  pressure;  LES/RANS  on 

Geometry  #2) 
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Figure  7:  Centerplane  reactant  mass-fraction 
contours 


4.3  Results  and  Discussion:  geometry 

Based  on  the  results  shown  above,  a  new  mesh  was  generated  that  maintains  good  resolution  in 
the  diffuser  section  while  resolving  all  injector  ports.  The  mesh  contains  ~229  M  cells 
partitioned  over  17676  blocks.  The  remaining  results  were  obtained  on  this  mesh.  Two  cases 
have  been  considered,  one  with  a  22.8  Torr  back  pressure  and  the  other  with  a  14.7  Ton- 
backpressure. 
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Figure  8:  Centerplane  Mach  number  contours  for  flow  with  22.8  Torr  back  pressure. 
Top:  RANS,  Middle:  time  averaged  LES/RANS,  Bottom:  Instantaneous  LES/RANS 

snapshot 
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To  initialize  the  LES/RANS  calculations,  the  two  cases  were  run  to  convergence  using  the 
Menter  SST  RANS  turbulence  closure  The  LES/RANS  simulation  with  22.8  Torr  back  pressure 
was  performed  using  a  timestep  of  1.5e-7  seconds  and  progressed  for  0.0135  seconds,  or 
approximately  nine  flow  residence  times  through  the  nozzle  and  diffuser  sections.  Time- 
averaged  statistics  were  collected  for  the  last  five  flow-through  times.  The  simulation  with  14.7 
Torr  back  pressure  was  performed  using  the  same  timestep  over  0.082  seconds  or  approximately 
eight  flow  residence  times,  with  statistics  gathered  over  the  last  three  flow  residence  times, 
which  is  not  yet  enough  for  statistical  convergence.  Both  cases  displayed  long-period 
oscillations  in  the  flow  resulting  in  continued  fluctuations  in  the  mass  flow  rates  through  the 
cavity  -  the  discrepancy  between  inflow  and  outflow  mass  fluxes  reached  2%  in  the  22.8  Torr 
case  and  as  high  as  7%  in  the  14.7  Torr  case. 
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Figure  9:  Centerplane  Ni  mass  fraction  contours  for  flow  with  22.8  Torr  back  pressure. 

Top:  RANS,  Middle:  time  averaged  LES/RANS,  Bottom:  Instantaneous  LES/RANS 

snapshot 

Figures  8-9  compare  the  centerplane  Mach  number  and  N2  mass  fraction  contours  from  the 
RANS,  instantaneous  LES/RANS,  and  time-averaged  LES/RANS  data  for  the  22.8  Torr  back 
pressure  case.  The  overall  structure  of  the  flow  for  this  case  can  be  seen:  a  series  of  bow  shocks 
due  to  the  injectors  help  disperse  the  injected  gas,  but  otherwise  the  flow  in  the  nozzle  remains 
supersonic  and  largely  coherent.  The  main  shock  train  forms  around  0.25  meters  downstream  of 
the  throat,  about  15  cm  upstream  of  the  supersonic  diffuser.  At  this  point,  the  core  flow 
decelerates  and  large-scale  turbulence  triggered  by  the  shock-boundary  layer  interaction  is 
observed  in  the  LES/RANS  snapshot.  A  feature  seen  in  the  LES/RANS  but  not  the  RANS  data 
is  the  large  separation  bubble  on  the  upper  wall  at  the  start  of  the  shock  train.  It  is  clear  from 
comparison  of  the  instantaneous  and  time-averaged  data  that  the  separated  region  is  a  transient 
feature,  but  its  residence  time  is  long  enough  to  affect  the  average. 
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Figure  10:  Centerplane  Mach  number  contours  for  flow  with  14.7  Torr  back  pressure. 
Top:  RANS,  Middle:  time  averaged  LES/RANS,  Bottom:  Instantaneous  LES/RANS 

snapshot 
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Figure  11:  Centerplane  N2  mass  fraction  contours  for  flow  with  14.7  Torr  back 
pressure.  Top:  RANS,  Middle:  time  averaged  LES/RANS,  Bottom:  Instantaneous 

LES/RANS  snapshot 

The  14.7  Torr  back  pressure  case  (Figures  10-1 1)  shows  similar  upstream  behavior  but  the  back 
pressure  does  not  propagate  through  the  supersonic  diffuser,  and  consequently  the  shock  train  is 
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more  gradual  and  comparatively  far  downstream.  In  fact  the  core  flow  remains  supersonic  as  it 
enters  the  subsonic  diffuser,  and  large-scale  turbulence  is  not  observed  until  that  point,  about  60 
cm  downstream  of  the  nozzle.  There  appears  to  be  a  slight  thickening  of  the  boundary  layer  on 
the  upper  wall  and  a  downward  deflection  in  the  flow,  but  it  is  not  as  severe  as  in  the  higher  back 
pressure  case  and  does  not  disrupt  the  flow  through  the  supersonic  diffuser  section. 


a.  22.8  Torr  back  pressure  b.  14.7  Torr  back  pressure 

Figure  12:  Comparison  of  upper  centerline  wall  pressures  for  different  values  of  the 

back  pressure 

A  comparison  of  the  computed  upper  wall  pressures  with  experimental  data  from  Noren,  et  al  (3) 
(Figure  12)  shows  that  in  both  cases  the  RANS  and  LES/RANS  simulation  are  accurately  able  to 
place  the  location  of  the  pressure  rise  associated  with  the  shock  train.  However,  in  the  22.8  Torr 
case  the  predicted  pressure  recovery  in  the  supersonic  diffuser  is  slower  than  it  is  in  the 
experiment,  for  both  the  RANS  and  LES/RANS  simulations.  The  slow  pressure  rise  is  also  seen, 
albeit  to  a  lesser  extent,  in  the  14.7  Torr  case.  This  suggests  an  under-prediction  of  the  boundary 
layer  thickening  due  perhaps  to  insufficient  mixing  with  the  free-stream,  which  will  be  examined 
later.  Overall,  the  predictions  from  the  RANS  and  LES/RANS  simulations  are  quite  similar, 
barring  a  pressure  dip  about  35cm  from  the  nozzle  in  the  LES/RANS  data  for  the  22.8  Torr  case. 
This  is  associated  with  the  transient  separation  bubble  noted  earlier  and  is  confirmed  by  the  spike 
in  the  pressure  variance  in  that  location.  The  standard  deviation  of  wall  pressure  in  both  cases 
signifies  high  levels  of  unsteadiness  within  the  shock  train  itself,  but  the  unsteadiness  is  much 
less  in  the  lasing  cavity  downstream  of  the  injection  region.  Figure  13  shows  streamwise 
velocity  snapshots  at  the  top  wall  in  the  22.8  Torr  case  at  points  where  the  shock  train  is  at  its 
maximum  and  minimum  upstream  locations.  The  blue  and  green  areas  show  where  the  flow  has 
detached.  While  there  is  variation  in  the  overall  size  of  the  regions  of  reversed  flow,  the  back 
end  of  the  shock  appears  to  remain  mostly  stationary,  with  reversed  flow  at  the  comers  reaching 
back  to  21cm  downstream  of  the  throat.  These  persistent  comer  regions  likely  provide  the 
pathway  for  the  upstream  pressure  propagation  into  the  lasing  cavity.  The  band  of  reversed  flow 


Approved  for  public  release;  distribution  unlimited. 

11 


Wall  pressure  standard  deviation  (torr) 


on  the  sidewalls  26em  downstream  of  the  throat  is  due  to  a  slight  expansion  in  the  ehannel  at  the 
start  of  the  diffuser. 
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Figure  13:  Maximum  (top)  and  minimum  (bottom)  extent  of  reversed  flow  on  the  upper 

wall  for  22.8  Torr  case 

A  key  question  in  the  understanding  of  this  flow  regime  is  the  level  of  interaetion  between  the 
freestream  and  boundary  layer  near  the  shock  train.  Figure  14  shows  snapshots  of  the  injected 
N2  isosurface  in  the  lasing  cavity  for  the  22.8  Torr  and  14.7  Torr  cases,  providing  a  visualization 
of  the  core  flow.  The  individual  injectant  plumes  appear  to  stay  coherent  far  downstream  of  the 
injectant  location.  Structural  changes  due  to  the  interaction  of  the  shock  train  with  the  plumes 
begin  to  occur  at  X  =  25  cm  for  the  22.8  Torr  case. 


a.  22.8  Torr  back  pressure 


b.  14.7  Torr  back  pressure 


Figure  14:  Snapshots  of  N2  mass  fraction  iso-surfaces  (0.065),  colored  by  pressure 
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As  the  shock  train  is  downstream  of  the  lasing  cavity  in  the  14.7  Torr  back-pressure  case,  no 
such  changes  are  evident.  Eddies  in  each  jet  are  likely  formed  by  Kelvin-Helmholtz  instability, 
which  could  be  enhanced  by  the  oblique  shocks  immediately  downstream  of  the  injectors. 
Cross-flow  contours  of  N2  mass  fraction  for  the  22.8  Torr  case  (Figure  15)  indicate  that  the 
actions  of  the  eddies  are  sufficient  to  almost  fully  mix  the  flow  by  the  time  it  reaches  the 
supersonic  diffuser.  Similar  results  (not  shown)  are  present  for  the  14.7  Torr  case.  Figure  15 
also  shows  the  instantaneous  Mach  number  in  the  cavity. 

Y 

Z 


Figure  15:  Instantaneous  N2  mass  fraction  and  Mach  number  profiles  in  lasing 

cavity 

Boundary  layer  growth  in  the  upstream  section  of  the  cavity  is  slow,  and  the  main  interaction 
between  the  boundary  layer  and  free  stream  appears  to  be  the  formation  of  persistent  vortices  at 
the  comers  which  entrain  some  slow-moving  gas  into  the  center  of  the  flow.  However, 
significant  boundary  layer  thickening  only  occurs  after  the  first  oblique  shock,  seen  in  the  fourth 
downstream  contour.  After  this,  the  supersonic  core  flow  rapidly  decelerates  and  the  visible  flow 
stmctures  de-cohere.  The  resolved  turbulence  kinetic  energy,  plotted  in  Figure  16  for  both  back¬ 
pressure  cases,  shows  that  there  is  little  to  no  turbulence  production  in  the  freestream,  which  is 
expected  given  the  low  Reynolds  numbers  involved  (the  present  case  has  a  Reynolds  number  of 
approximately  40000  based  on  the  nozzle  length).  Turbulence  energy  is  only  present  close  to  the 
wall  and  is  not  significantly  mixed  by  the  free-stream  eddies,  even  at  the  comers.  Interaction 
with  the  shock  train  causes  increased  turbulence  production,  especially  within  the  separated  flow 
region  at  the  top  wall,  but  this  energy  dissipates  farther  into  the  shock  train.  In  the  farthest 
downstream  contour  we  can  see  that  the  concentration  of  turbulence  energy  remains  very  close  to 
the  wall  with  a  limited  amount  mixed  into  the  core  flow  eddies.  This  energy  concentration 
causes  a  continued  thinning  of  the  boundary  layer,  explaining  the  overly  slow  predicted  wall 
pressure  rise  in  the  supersonic  diffuser. 


b.  Macb  number 
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a.  N2  mass  fraction  b.  Mach  number 

Figure  16:  Resolved-scale  turbulence  kinetic  energy  in  the  lasing  cavity 


Overall,  the  limited  interaction  between  the  core  flow  and  boundary  layers  means  that  accurate 
modeling  of  the  flow  will  be  highly  sensitive  to  the  turbulence  model’s  ability  to  capture 
boundary  layer  behavior  and  shock-boundary  interaction.  Additionally,  because  the  bulk  of  the 
turbulence  production  occurs  near  the  wall  where  the  LES/RANS  model  transitions  to  RANS, 
this  explains  why  the  wall  pressure  prediction  between  the  pure  RANS  and  hybrid  models  is  so 
similar. 

5.0  CONCLUSIONS 

Several  hybrid  LES/RANS  and  RANS  simulations  have  been  conducted  for  model  COIL 
laser  flow  fields  on  meshes  with  different  levels  of  isotropic  resolution  and  under  different  back 
pressures.  Adequate  resolution  of  outer-layer  turbulent  eddies  is  of  primary  importance  in 
ensuring  reasonable  success  of  the  LES/RANS  modeling  strategy.  Estimates  of  outer-layer 
turbulent  length  scales  obtained  from  pre-cursor  RANS  simulations  have  been  used  to  design 
large  meshes  that  provide  good  resolution  for  capturing  the  transition  to  turbulence  induced  by 
the  shock  train.  Wall-pressure  predictions  obtained  by  both  LES  and  LES/RANS  indicate 
reasonable  agreement  with  experimental  data,  with  the  most  notable  discrepancy  being  the 
inability  of  either  model  to  capture  the  rate  of  pressure  increase  experienced  as  the  flow  passes 
through  the  leading  edge  of  the  shock  train.  The  LES/RANS  model  is  capable  of  predicting 
shock-induced  transition  to  turbulence  as  well  as  unsteadiness  due  to  instability  growth  as  the 
reactant  jets  enter  and  mix  within  the  lasing  cavity.  Reactant  mixing  is  due  to  laminar- flow 
effects  and  is  performed  efficiently  by  the  jet-array  configuration.  Local  separation  regions 
induced  by  the  shock  train  do  not  propagate  very  far  into  the  mixing  region,  except  possibly  at 
the  highest  back  pressure  of  26  Torr  (which  was  simulated  on  a  mesh  that  did  not  contain  the 
injection  array).  The  effects  of  shock-train  propagation  on  reactant  mixing  appear  to  be  minor  at 
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lower  back  pressures  (14.7  Torr  and  22.8  Torr),  indicating  that  the  current  COIL  diffuser  design 
might  provide  adequate  isolation. 
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LIST  OF  SYMBOLS,  ABBREVIATIONS  AND  ACRYNOMS 


LES/RANS 

Large-eddy  Simulation/  Reynolds-averaged  Navier-Stokes 

COIL 

Chemical  Oxygen  Iodine  Laser 

RANS 

Reynolds- Averaged  Navier-Stokes 

NCSU 

North  Carolina  State  University 

ILU 

Incomplete  Lower-Upper 

LDFSS 

Low  Diffusion  Flux- Splitting  Scheme 

LD-PPM 

Low-Dissipation  Piecewise  Parabolic  Method 

PPM 

Piecewise  Parabolic  Method 

SST  RANS 

Shear  Stress  Transport  RANS 
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